'''
Created on 20/01/2013

@author: Jorge
'''
import numpy as np

class OptimizateProbabilityParam(object):
    '''
    classdocs
    '''


    def __init__(self, model, X, Y):
        '''
        Constructor
        '''
        self.model = model
        #self.X = X
        
        self.f = self.model.predict_proba(X)[:,1]
        N_pos = len([y for y in Y if y==1])
        N_neg = len(Y)-N_pos
        t = []
        for item in Y:
            if item==1:
                t.append( (N_pos+1.0)/(N_pos+2.0))
            else:
                t.append( 1.0/(N_neg+2.0))
        self.t = np.array(t)  #target probability


    #params = A and B
    def minimaze(self, param):
        A, B = self.get_param(param)
        fx = self.f
        pk = self.sigmoid(A,fx,B)
        #tk = self.t
        #pk = 1.0/( 1+ np.exp(A[j]*self.f[k]+B[j]) )
        #suma = np.dot(self.t, np.log(pk)) + np.dot( (1-self.t), np.log(1-pk) )
        #print 'pk= ',pk
        suma=0
        print 'A,B=', A, B
        for k in range(len(self.t)):
            if self.t[k]==0:
                suma+= np.sum( (1-self.t[k])*np.log(1-pk[k]) )
            elif self.t[k]==1:
                suma+= np.sum( self.t[k]*np.log(pk[k]) )
            else:
                suma+= np.sum( self.t[k]*np.log(pk[k])+(1-self.t[k])*np.log(1-pk[k]) )
        suma = suma*-1
        if suma<1.22124532709e-13:
            suma=0
        print suma
        return suma
    
    #devolver vector de 24 elementos
    def gradient(self, param):
        A, B = self.get_param(param)
        grad= np.zeros([2])
        grad[0]= self.fprime_A(A, B)
        grad[1]= self.fprime_B(A, B)
        
        return grad
    
    def fprime_A(self, A, B):
        fx = self.f
        pk = self.sigmoid(A,fx,B)
        tk = self.t
        return np.sum( (pk-tk)*fx )
    
    def fprime_B(self, A, B):
        fx = self.f
        pk = self.sigmoid(A,fx,B)
        tk = self.t
        return np.sum( pk-tk )
            
    def sigmoid(self, A, f, B):
        return 1.0/( 1+ np.exp(A*f+B) )
    
    def hessian_matriz(self, param):
        A, B = self.get_param(param)
        matriz= np.zeros([2,2])
        matriz[0,0] = self.deltaA_deltaA(A, B)
        matriz[0,1] = self.deltaA_deltaB(A, B)
        matriz[1,0] = self.deltaB_deltaA(A, B)
        matriz[1,1] = self.deltaB_deltaB(A, B)
        
        return matriz
    
    def deltaA_deltaA(self, A, B):
        fx = self.f
        pk = self.sigmoid(A,fx,B)
        tk = self.t
        return np.sum( np.power(fx,2)*pk*(pk-1))
    
    def deltaA_deltaB(self, A,B):
        fx = self.f
        pk = self.sigmoid(A,fx,B)
        tk = self.t
        return np.sum( fx*pk*(pk-1) )
    
    def deltaB_deltaB(self, A,B):
        fx = self.f
        pk = self.sigmoid(A,fx,B)
        tk = self.t
        return np.sum( pk*(pk-1) )
    
    def deltaB_deltaA(self, A,B):
        fx = self.f
        pk = self.sigmoid(A,fx,B)
        tk = self.t
        return np.sum( fx*pk*(pk-1) )
    
    
    def get_param(self, param):
        return param[0], param[1]
        
    